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LONG-TERM  GOAL 

The  primary  focus  of  this  research  is  to  use  large-eddy  simulation  (LES)  and  large-wave  simulation 
(LWS)  to  obtain  improved  physical  understanding  of  wind-wave-ocean  interactions,  based  on  which 
we  aim  to  develop  effective  models  of  wind  input  and  whitecapping  dissipation  for  phase-resolving, 
nonlinear  wave-field  simulation  at  large  scales.  Our  ultimate  goal  is  to  establish  a  numerical  capability 
for  predicting  deterministically  large-scale  nonlinear  wave-field  in  real  marine  environments  with  the 
presence  of  significant  wind  and  wave  breaking  effects. 

OBJECTIVES 

The  scientific  and  technical  objectives  of  this  research  are  to: 

•  develop  advanced  LES  and  LWS  numerical  capabilities  for  wind- wave-ocean  interactions  with 
physics-based  subgrid-scale  (SGS)  models;  use  high-performance  LES/LWS  as  a  powerful 
research  tool  to  obtain  an  improved  understanding  of  the  flow  structure  in  the  atmosphere-ocean 
wave  boundary  layer 

•  develop  effective  models  for  wind  input  and  the  associated  whitecapping  dissipation  in  a  direct 
phase-resolving  context,  which  can  be  readily  incorporated  into  the  deterministic  numerical  tool  of 
the  Simulation  of  Nonlinear  Ocean  Wave-field  (SNOW) 

•  understand  effects  of  multi-scale  physics  and  environmental  uncertainties  upon  wave  deterministic 
propagation,  and  to  effectively  model  these  effects;  validate  the  direct  modeling  and  simulation 
approach,  and  perform  direct  comparison  with  existing  theories  and  field  measurements 

APPROACH 

We  use  a  systematic,  multiscale  approach  to  investigate  and  to  model  effects  of  wind  input  and 
whitecapping  dissipation  on  wave-field  evolution.  This  includes:  (1)  use  LES  and  LWS  to  obtain 
improved  physical  understanding  of  wind-wave-ocean  interactions  at  small  scales  (0(1  ~10)  significant 
gravity  wave  lengths);  (2)  based  on  the  LES/LWS  results,  develop  advanced  wind  input  and 
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whitecapping  models  in  a  direct  physical  context  in  terms  of  surface  pressure  distribution  and  flow 
field  filtering,  respectively;  and  (3)  incorporate  the  models  into  SNOW  simulation  to  investigate  local 
effects  of  wind  forcing  and  whitecapping  on  large-scale  wave-field  evolution.  Because  the  physics  are 
being  investigated  and  modeled  in  a  direct,  phase-resolving  context,  we  expect  that  the  models 
developed  in  this  study  are  more  likely  to  succeed  than  traditional  phase-averaged  parameterizations. 

The  numerical  study  of  wind-ocean-wave  interactions  are  at  two  fronts:  viscous  flow  simulation  for 
turbulence-wave  interactions  at  small  scales,  and  potential  flow  based  wave  simulation  at  large  scales. 
For  viscous  air-wave-water  simulations,  we  use  the  approach  of  large-eddy  simulation  based  on 
filtered  Navier-Stokes  equations,  in  which  large  turbulence  eddies  are  computed  explicitly  with  effects 
of  small  eddies  being  represented  by  SGS  models.  A  similar  LWS  approach  is  used  for  wave- 
turbulence  interactions,  with  large  wave  components  being  simulated  directly  and  small  wave  effects 
modeled. 

One  of  the  major  issues  with  simulating  coupled  air- wave- water  turbulent  flows  is  the  presence  of  a 
deformable,  time-evolving  free  surface,  which  makes  the  air  and  water  computational  domains 
irregular.  The  location  and  geometry  of  the  free  surface  are  unknown  beforehand,  and  they  are  part  of 
the  solution  to  be  solved  for.  To  overcome  this  difficulty,  we  have  developed  a  suite  of 
complementary  computational  methods,  which  include  a  boundary  interface  tracking  method  based  on 
boundary-fitted  grids  for  moderate  wave  slope  and  an  Eulerian  interface  capturing  method  based  on  a 
level  set  approach  for  steep/breaking  waves. 

The  wind  input  and  whitecapping  models  developed  from  LES/LWS  study  will  be  incorporated  to 
potential  flow  wave  computation  using  SNOW  developed  at  MIT.  The  SNOW  uses  a  high-order 
spectral  (HOS)  method,  which  is  a  pseudo- spectral  method  developed  based  on  the  Zakharov  equation 
and  mode-coupling  idea.  The  wind  input  to  the  wave-field  is  modeled  by  surface  distribution  of 
pressure,  while  wave  breaking  dissipation  is  represented  by  advanced  filtering  treatment  in  physical 
and  spectral  spaces. 

The  multiscale  modeling  of  ocean  wave-fields  will  be  one  of  the  foci  of  this  research.  At  local  scales, 
we  apply  LES  and  LWS,  with  which  the  large  eddies  and  waves  are  resolved  explicitly  and  structures 
smaller  than  the  computational  grid  are  modeled.  The  wind  input  and  models  developed  based  on 
local-scale  LES/LWS  are  then  implemented  in  SNOW  for  large-scale  simulations.  In  SNOW,  the 
nonlinear  wave  interactions,  the  wind  input,  and  the  whitecapping  dissipation  are  all  represented  in  a 
direct  physical  context.  This  approach  provides  a  powerful  tool  and  a  unique  opportunity  to 
investigate  the  effects  of  local  processes  on  large-scale  wave-field  evolution  in  a  physical-based, 
deterministic  way. 

Large-scale  high-performance  computation  on  parallel  computers  is  used  to  meet  the  computational 
challenges  in  the  turbulence  and  wave  simulations.  Massage  passing  interface  (MPI)  based  on  domain 
decomposition  is  used  for  parallelization. 

WORK  COMPLETED 

During  the  fiscal  year  of  2009,  substantial  progresses  have  been  made,  which  include: 
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•  We  have  validated  our  simulation  results  through  extensive  comparison  with  measurement  and 
simulation  data  in  the  literature. 

•  We  have  quantified  and  characterized  wind  input  to  nonlinear  irregular  wavefield  that  evolves 
dynamically  under  wind  action  in  the  simulation. 

•  We  have  analyzed  the  pressure  variation  near  the  wave  surface,  the  result  of  which  is  helpful 
for  field  measurement  of  wind  input  to  waves. 

•  We  have  extended  the  simulation  capability  to  turbulence  transport  of  scalars  in  marine 
atmospheric  boundary  layer. 

•  We  have  performed  direct  simulation  of  steep  and  breaking  wind-waves  and  airflow  separation, 
and  quantified  the  associated  wind-wave  interaction  processes. 

•  We  have  extended  the  simulation  capability  by  developing  a  hybrid  simulation  approach  using 
state-of-the-art  Eulerian  and  Lagrangian  methods  for  free-surface  flow  simulation. 

RESULTS 

In  previous  years  of  this  project,  we  developed  novel  numerical  methods  for  wind  and  wave  motions; 
from  simulations,  the  fundamental  dynamics  of  wind-wave  interactions  were  investigated.  In  FY  2009, 
we  continued  to  investigate  wind  input  to  wave  evolution  and  performed  extensive  analysis  to  connect 
our  results  to  other  studies.  Figure  1  shows  the  dependence  of  wave-induced  form  drag  on  wave  age 
c/U  a/2.  When  the  wave  age  is  small,  the  form  drag  has  a  positive  value,  which  results  in  growth  of  the 
water  wave.  As  the  wave  age  increases,  the  form  drag  increases  first  and  reaches  its  peak  value,  and 
then  decreases  as  the  wave  age  further  increases.  The  form  drag  vanishes  around  c/U  a/2  =  1,  and 
becomes  negative  for  large  wave  age,  which  results  in  wave  attenuation.  The  simulation  result 
obtained  in  this  study  is  consistent  with  other  simulation  results  in  the  literature,  which  are  plotted  in 
figure  1  for  comparison. 

Figure  2  shows  the  dependence  of  wave  growth  rate  parameter  0  on  the  wave  steepness  ak  for  the 
slow  wave  case.  As  wave  steepness  increases,  the  growth  rate  parameter  decreases.  The  measurement 
data  compiled  by  Peirson  &  Garcia  (2008)  is  also  plotted  for  comparison.  It  is  found  that  the  result  of 
the  present  study  agrees  well  with  the  measurements  in  the  literature. 

Under  the  wind  pressure  forcing,  the  temporal  growth  rate  of  a  water  wave  can  be  parameterized  as 
Y  =  S  a/2  (Ua/2/c-1)|  U  a/2/c-1|,  where  S  x/2  is  the  sheltering  coefficient.  Based  on  laboratory 
measurement,  Donelan  (1999)  obtained  S  a/2  =  0.28.  Donelan  et  al.  (2006)  further  improved  this 
parameterization  by  using  wave-following  field  measurement  data  and  showed  that  S  a/2  =  0.17.  Figure 
3  shows  the  comparison  between  our  simulation  result  and  the  parameterizations.  It  is  found  that  our 
simulation  result  agrees  very  well  with  the  measurement  and  supports  their  improvement  in  the 
parameterization. 

To  investigate  the  interaction  between  wind  turbulence  and  broadband  wavefield,  we  performed  a 
series  of  turbulence-wave  coupled  simulations.  We  constructed  the  initial  wavefield  based  on  the 
JONSWAP  spectrum.  Three  cases  with  different  wave  ages  based  on  the  phase  velocity  at  the  peak 
wavenumber  km,  cm/u*=5,  12.3,  and  16,  are  considered.  Figure  4  shows  the  wave  growth  rate 
parameter  p  and  the  fractional  rate  of  energy  input  y  as  a  function  of  wavenumber  k  for  these  three 
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cases.  To  help  understand  their  behavior,  we  also  plot  the  values  of  c/u*  at  different  k  for  the  three 
cases.  Note  that  according  to  the  wave  dispersion  relationship,  c/u*  decreases  as  k  increases  and  has 
different  ranges  for  the  three  cases.  For  the  case  of  cm/u*=5,  P  reaches  its  peak  around  k  =  km  and 
decreases  as  k  further  increases.  For  the  cases  of  cm/u*=12.3  and  16,  P  first  increases  with  k  and 
reaches  its  peak  around  k  =  2km~3km,  and  then  decreases  as  k  further  increases.  For  all  cases,  p 
reduces  significantly  at  large  k  (i.e.,  short  waves  in  the  spectrum);  at  these  small  scales,  the  pressure 
induced  by  the  short  waves  is  relatively  small  compared  to  the  pressure  fluctuation  in  the  wind 
turbulence.  On  the  other  hand,  y,  which  measures  the  fractional  rate  of  energy  input  to  wave,  does  not 
reduce  as  rapidly  compared  to  p  (y  even  increases  with  k  for  the  case  of  cm/u*=5).  In  reality,  local 
small  waves  may  grow  rapidly  to  break,  serving  as  an  important  vehicle  for  atmosphere-ocean 
momentum  and  energy  transfer. 

We  also  extended  our  study  to  turbulence  transport  of  scalars  in  the  wind- wave  field.  Figure  5  shows 
the  phase-averaged  contours  of  the  scalar  concentration  fluctuation  root-mean- square  (rms)  value 
above  water  waves  with  wave  ages  c/u*=2  and  14.  For  the  case  of  c/u*=2,  the  scalar  fluctuation  is 
high  in  the  strong  shear  layer  that  starts  from  the  wave  crest  to  above  the  wave  trough.  For  the  case  of 
c/u*=14,  the  high  intensity  region  of  scalar  fluctuation  is  located  above  the  windward  face  of  the  wave 
crest.  This  dependence  of  scalar  fluctuations  on  the  wave  age  is  found  to  be  similar  to  that  of  the 
streamwise  velocity  fluctuation. 

For  accurate  estimation  of  the  surface  pressure  in  field  measurement,  it  is  necessary  to  put  the  pressure 
probe  close  to  the  wave  surface.  In  figure  6,  we  compare  the  vertical  pressure  decay  profiles  of  our 
slow  wave  cases  to  the  Kudyavtsev  et  al.  (2001)  pressure  profile  shown  by  Donelan  et  al.  (2006)  with 
their  record  LG9.  Two  exponential  profiles,  exp(-kz)  and  exp(-0.55kz)  are  also  plotted  for  comparison. 
The  result  shows  that  the  pressure  decay  rate  decreases  as  the  wave  steepness  increases.  When  the 
measurement  is  taken  at  kz=0.2,  the  surface  pressure  values  from  our  simulations  collapse,  and  the  two 
exponential  extrapolations  provide  values  with  about  5%  error;  when  the  measurement  is  taken  at  kz=l 
instead,  the  surface  pressure  values  from  different  cases  of  our  simulations  diverge,  and  the 
extrapolations  can  have  errors  as  large  as  20%.  This  comparison  confirms  the  advantage  of  wave¬ 
following  pressure  probe  measurement. 

Under  high  wind  conditions,  the  wind-wave  becomes  steep  and  wave  breaking  may  happen,  which 
results  in  airflow  separation  and  wave  energy  dissipation.  Figure  7  shows  the  flow  structures  when 
instantaneous  airflow  separation  happens  over  steep/breaking  wind-wave.  The  enhanced  airflow  shear 
layer  at  the  edge  of  the  separation  zone  is  clearly  indicated  by  the  high  intensity  of  spanwise  vorticity. 
The  violent  turbulence  motion  around  the  separation  zone  intensifies  the  turbulence  mixing  and 
momentum  transport.  Our  analysis  also  shows  that  the  strong  reattachment  of  the  airflow  significantly 
increases  the  surface  pressure  and  shear  stress  on  the  windward  face  of  the  preceding  wave  crest,  and 
therefore  enhances  the  form  drag  and  tangential  stress  on  the  wave  surface. 

When  wave  breaking  happens,  the  surface  process  becomes  violent  and  full  of  small  scale  structures. 

In  order  to  study  the  details  of  wave  breaking  as  well  as  the  associated  airflow  separation  and  spray 
droplet  dynamics,  we  developed  a  novel  hybrid  simulation  tool  by  coupling  the  strength  of  several 
multiphase  fluid  simulation  methods.  We  use  the  coupled  level-set  and  volume-of-fluid  (CLSVOF) 
method  to  simulate  the  wind- wave  interaction  and  the  development  of  breaking  wave  at  the  full  scale. 
Once  wave  breaking  happens,  locally  a  high  resolution  smoothed  particle  hydrodynamics  (SPH) 
simulation  is  coupled  with  CLSVOF  to  resolve  and  track  the  violent  motion  of  the  fluid  particles. 

Figure  8  shows  an  example  of  this  CLSVOF-SPH  coupled  simulation. 


4 


IMPACT/APPLICATION 


This  project  aims  at  a  basic  scientific  understanding  of  the  air-sea- wave  interaction  physics  and 
numerical  capability  development  for  ocean  wave-field  prediction  pertinent  to  Navy  applications.  It 
addresses  a  critical  need  of  the  Navy  to  bridge  the  gap  between  the  modeling  of  small-scale  air-sea- 
wave  interaction  physics  and  the  prediction  of  ocean  waves  at  regional  scales.  The  proposed  research 
will  provide  the  Navy  with  a  new  powerful  tool  to  predict  deterministically  nonlinear,  large  wavefield 
with  finely-resolved  temporal  and  spatial  details.  The  new  phase-resolved,  deterministic  tool  is 
fundamentally  distinct  from  existing  phase-averaged,  statistical  wave  modeling  tools  such  as  WAM 
and  SWAN,  with  the  potential  of  being  able  to  make  more  accurate  prediction  because  of  its  direct, 
physics-based  approach.  Furthermore,  the  results  of  this  work  will  also  be  useful  for  the  comparison 
and  calibration  of  field  measurements  and  for  obtaining  physical  insights  to  improve  existing  phase- 
averaged  wave  prediction  models. 

RELATED  PROJECTS 

This  work  complements  a  number  of  on-going  ONR  projects.  In  particular,  it  is  closely  related  to  the 
development  of  Simulation  of  Nonlinear  Ocean  Wave-field  (SNOW)  by  Professor  Dick  Yue’s  research 
group  at  MIT.  The  wind  input  and  wave  breaking  dissipation  modeling  in  this  project  is  to  be 
incorporated  into  SNOW,  and  together  we  will  improve  the  SNOW  capability  for  it  to  become  a  next 
generation  of  wave  model  capable  of  predicting  nonlinear  wave  evolution  subject  to  winds  and 
whitecapping.  Such  numerical  tools  will  be  useful  for  high-resolution  wave-field  study. 
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Figure  1.  Wave-induced  form  drag  Fp  as  a  function  of  c/  U  xn  for  water  wave  with  steepness 
ak=0.1.  A,  current  DNS  results  O,  DNS  results  of  Sullivan  et  al.  (2000)  O,  DNS  results  of 
Kihara  et  al.  (2007)  V,  results  of  Mastenbroek  et  al.  (1996)  with  LRR  closure.  Results  of  Li  et  al. 
(2000)  with  different  turbulent  closure:  —  #  — ,  E-  xz;  —  ** — >  q2l; - ,  LRR. 
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Figure  2.  Wave  growth  rate  parameter  &  for  slow  water  waves  plotted  as  a  function  of  wave 
steepness  ak.  Measurement  data  compiled  by  Peirson  &  Garcia  (2008)  are  denoted  by  solid  symbols 
with  error  bars:  A,  Bole  (1967);  U,  Banner  (1990);  ▼,  Mastenbroek  et  al.  (1996);  ♦,  Banner 

&  Peirson  (1998);  and  #,  Peirson  &  Garcia  (2008). 

Results  of  the  present  study  are  denoted  by  A . 
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Figure  3.  Comparison  of  temporal  growth  rate  Y  with  Donelan  et  al.  (2006):  — ,  their 
parameterization  Y=0.1 7(U x/2/c-l) \  U \n/c-l\;  V,  their  hindcast  estimates  based  on  the  data  of 

Donelan  (1999); - ,  parameterization  of  Donelan  (1999)  Y=0.28(U  x/2/c-l)\U  x/2/c-l\.  Results  of 

current  simulation  are  denoted  by  open  symbols:  □,  Airy  wave  with  various  steepness;  A,  Airy 
wave  with  steepness  ak=0.25;  O,  Stokes  wave  with  ak=0.25. 
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Figure  4.  (a)  Wave  growth  rate  parameter  p,  and  (b)  fractional  rate  of  energy  input  y,for  three 
broadband  wave  cases:  A,  cm/u*=5;  □,  cm/u*=12.3;  and  V,  cm/u *=1 6. 
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Figure  5.  Contours  of  phase-averaged  scalar  fluctuation  rms  above  water  waves  with 
steepness  ak=0.25  and  wave  ages:  (a)  c/u*=2;  and  (b)  c/u*=14.  The  surface  wave  propagates 

in  the  +x-direction. 


Figure  6.  Comparison  of pressure  profiles  with  exponential  decays:  solid  line,  the  Kudryavtsev  et  al. 
(2001)  profile  estimated  by  Donelan  et  al.  (2006)  with  steepness  ak=0.05  and  wave  age  c/u*=3.5 
(their  record  LG9);  dashed  line,  present  result  for  water  wave  with  (ak,  c/u*)  =  (0.1,  2);  dash-dot 
line,  the  present  result  for  water  wave  with  (ak,  c/u*)  =  (0.25,  2);  dash-dot-dot  line,  exponentially 
decay  with  exp(-kz);  dotted  line,  exponentially  decay  with  exp(-0.55kz).  The  profiles  are  normalized 
by  the  pressure  amplitude  \p\  at  (a)  kz—0.2  and  (b)  kz=l,  respectively  (i.e.,  assuming  that  the 
measurements  are  conducted  at  the  corresponding  heights). 
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Figure  7.  Instantaneous  airflow  separation  over  a  steep  wave  prior  to  breaking  under  high  wind 
condition  with  two  different  reattachment  condition  on  the  windward  face  of  the  preceding  wave 
crest:  (a)  weak  reattachment,  and  (b)  strong  reattachment.  The  instantaneous  three-dimensional 
velocity  vectors  are  projected  and  plotted  on  a  (x;  z)-plane  crossing  the  center  of  the  separation  zone. 

The  flow  fields  are  plotted  in  a  frame  following  the  wave  crest.  The  strong  shear  layers  at  the 
edge  of  the  separation  bubble  and  near  the  reattachment  point  are  indicated  by  the  high  intensity 

of  spanwise  vorticity. 
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Full  scale  Eulerian  simulation  using  coupled 
level-set  and  volume-of-fluid  method 


Local  scale  high  resolution 
Lagrangian  simulation  using 
smoothed  particle  hydrodynamics 


Figure  8.  CLSVOF  and  SPH  coupled  simulation. 
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